From Risk to Resilience: Ecological Exploration of Psychopathology and Flourishing in Neurodivergent Children Facing Adversity
This study extends analyses conducted by Bethell et al. (2023), who explored the role of family resilience and connection in promoting flourishing among U.S. children facing adversity. Our study applies a similar framework but focuses on neurodivergent children, examining how adverse childhood experiences (ACEs) influence mental health outcomes and how flourishing behaviors—such as curiosity, emotional control, and task persistence—may buffer these effects. This study was conducted by the Slopen Laboratory at the Harvard T.H. Chan School of Public Health.
Author: Adrián Medina
Date: December 5, 2024
Project Overview
This study investigates the relationship between adverse childhood experiences (ACEs) and mental health outcomes in children with neurodevelopmental disorders (NDDs), such as Autism Spectrum Disorder (ASD) and ADHD. The study specifically examines the moderating role of child flourishing behaviors—curiosity in learning, emotional control, and task persistence—on the effects of ACEs on mental health outcomes.
Study Design:
The study utilizes data from the National Survey of Children’s Health (NSCH), focusing on children aged 6-17 with reported neurodevelopmental issues (N = 54,528, MAge = 12.1). This cross-sectional analysis assesses the exposure to ACEs and its impact on three key mental health outcomes: anxiety, depression, and behavioral problems. Interaction terms in logistic regression models are used to assess whether child flourishing behaviors buffer the effects of ACEs on these outcomes.
Data Collection:
Key variables collected include:
Adverse Childhood Experiences (ACEs): Derived from parent-reported experiences such as exposure to violence, household dysfunction, and discrimination.
Mental Health Outcomes: Parent-reported diagnoses of anxiety, depression, and behavioral problems, confirmed by health care providers or educators.
Child Flourishing Metrics: Based on parent-reported behaviors such as curiosity in learning, emotional control, and task persistence, which were used to develop a Child Flourishing Index (CFI).
Analytic Approach:
Logistic regression models were employed to evaluate the association between ACEs and mental health outcomes, including interaction terms to investigate the moderating role of child flourishing. The models incorporated a comprehensive set of covariates, including child age, sex, race/ethnicity, socioeconomic status, caregiver mental health, neighborhood conditions, and healthcare access metrics.
Goals:
The study aims to:
Examine the dose-response relationship between ACE exposure and mental health outcomes among neurodivergent children.
Assess the protective role of child flourishing behaviors, with the goal of identifying resilience-promoting factors that can inform interventions for this vulnerable population.
Tip
For detailed information on the NSCH-derived variables, refer to the NSCH Data Dictionary.
Set up the R environment by configuring the CRAN repository, installing essential packages, setting the base path, and loading the primary dataset into the local environment.
Expand Code
# Set CRAN repository for consistent package installationoptions(repos =c(CRAN ="http://cran.rstudio.com/"))# Load or install the pacman package for efficient package managementif (!require(pacman)) install.packages("pacman")# Use p_load to load or install the required packagesp_load( dplyr, # Data manipulation tidyr, # Data reshaping (e.g., pivot_longer) readr, # Reading data (CSV files) knitr, # Creating tables for reports kableExtra, # Enhancing `knitr::kable()` tables ggplot2, # Data visualization sjlabelled, # Labelled data utility functions sjmisc, # Additional data manipulation functions survey, # Complex survey analysis lmtest, # Testing models (Breusch-Pagan test) and robust standard errors sandwich, # Calculating robust standard errors gtsummary, # Creating regression summary tables ggeffects, # Creating prediction plots gt, # Adjusting frequency/contingency tables psych, # Factor analysis crosstable, # Cross-Tabulations ggpubr, # Summary statistics + visualizations tidyverse, # Data visualization supplement ggdist # Distribution visualizations )# Specify the 'base_path' where you can find your data files, ASSUMING they're in the same directory, & set it as WDbase_path <-"~/GitHub/Adversity-Neurodiversity/data_files"setwd(base_path)# Load primary data fileneurodivergent_data <-read_csv("~/Downloads/nsch_2023e_topical_SAS/neurodivergent_data_2016to2023.csv")
1
See details under the Data Extraction & Filtering (Archived Code) callout below!
Data Extraction & Filtering (Archived Code)
The file being used for these analyses is a subset of a “master” file, Data2_2016to2022.csv, containing the compiled contents of the data releases from the National Survey on Child Health’s 2016-2022 cycles. The master file is >500 MB (0.5 GB) in size, making it cumbersome to process in real-time. For efficiency and ease of analysis, I am using a smaller subset that contains only the relevant data. This subset was generated using the data extraction and filtering process detailed in the archived code below, in an effort to maintain transparency and ensure reproducibility of the workflow.
Additionally, data from the NSCH 2023 was also later merged into the subset file being used for this analyses. The code for the extraction/wrangling process is shown below.
Expand Code
# Define the path to the datasetdata_path <-"~/Downloads/Data2_2016to2022.csv"Data2_2016to2022 <-read.csv(data_path)# Select specific variables, create 'neurodivergent' variable, & filter data based on age and neurodivergent statusneurodivergent_data <- Data2_2016to2022 %>%select(year, fpl, fpl_mean, sc_age_years, sc_hispanic_r, sc_race_r, sc_sex, higrade_tvis, k2q35a, k2q35c, k2q38a, k2q38c, k2q30a, k2q30c, k2q36a, k2q36c, k2q60a, k2q60c, k2q37a, k2q37c, downsyn, downsyn_desc, k2q31a, k2q31c, k2q61a, cerpals_desc, k2q33a, k2q33b, k2q33c, k2q32a, k2q32b, k2q32c, k2q34a, k2q34b, k2q34c, ace1, ace3, ace4, ace5, ace6, ace7, ace8, ace9, ace10, ace11, ace12, k6q71_r, k7q85_r, k7q84_r, talkabout, wktosolve, strengths, hopeful, k8q21, k8q30, fwc, k2q42a, k2q42c, a1_menthealth, a2_menthealth, k2q01, k4q27, k4q28x01, k4q28x04, appointment, available, issuecost, notelig, notopen, transportcc, family, family_r, k10q11, k10q12, k10q13, k10q14, k10q20, k10q22, k10q23, k10q40_r, k10q41_r, k10q30, k10q31, k9q96, sc_cshcn, k3q22, k3q20, instype, treatneed, k4q22_r, menbevcov) %>%# Creating 'neurodivergent' variable as NSCH does not explicitly inquire about neurodivergent status.# Using a list of diagnoses provided by NSCH to define neurodivergent individuals.mutate(neurodivergent =if_else(k2q35a ==1| k2q38a ==1| k2q36a ==1| k2q60a ==1| k2q37a ==1| k2q30a ==1| downsyn ==1| k2q31a ==1| k2q61a ==1| k2q42a ==1, 1, 0)) %>%# Filter data to include only individuals aged 6 or above and identified as neurodivergentfilter(sc_age_years >=6, neurodivergent ==1)# Resulting filtered data to be used for further analysiswrite.csv(neurodivergent_data, "neurodivergent_data_2016to2022.csv")
Expand Code
# Define the path to the datasetdata_path <-"~/Downloads/nsch_2023e_topical_SAS/nsch_2023e_topical.csv"nsch_data <-read.csv(data_path)# Create a new variable 'fpl_mean' representing the average of specific columnsnsch_data <- nsch_data %>%mutate(fpl_mean =rowMeans(select(., fpl_i1, fpl_i2, fpl_i3, fpl_i4, fpl_i5, fpl_i6), na.rm =TRUE))# Select specific variables, create 'neurodivergent' variable, & filter data based on age and neurodivergent status (manually added in blank columns fpl and downsyn_desc)# `menbevcov` was not collected for 2023neurodivergent_data_23 <- nsch_data %>%select(year, fpl, fpl_mean, sc_age_years, sc_hispanic_r, sc_race_r, sc_sex, higrade_tvis, k2q35a, k2q35c, k2q38a, k2q38c, k2q30a, k2q30c, k2q36a, k2q36c, k2q60a, k2q60c, k2q37a, k2q37c, downsyn, downsyn_desc, k2q31a, k2q31c, k2q61a, cerpals_desc, k2q33a, k2q33b, k2q33c, k2q32a, k2q32b, k2q32c, k2q34a, k2q34b, k2q34c, ace1, ace3, ace4, ace5, ace6, ace7, ace8, ace9, ace10, ace11, ace12, k6q71_r, k7q85_r, k7q84_r, talkabout, wktosolve, strengths, hopeful, k8q21, k8q30, fwc, k2q42a, k2q42c, a1_menthealth, a2_menthealth, k2q01, k4q27, k4q28x01, k4q28x04, appointment, available, issuecost, notelig, notopen, transportcc, family_r, k10q11, k10q12, k10q13, k10q14, k10q20, k10q22, k10q23, k10q40_r, k10q41_r, k10q30, k10q31, k9q96, sc_cshcn, k3q22, k3q20, instype, treatneed, k4q22_r) %>%# Creating 'neurodivergent' variable as NSCH does not explicitly inquire about neurodivergent status.# Using a list of diagnoses provided by NSCH to define neurodivergent individuals.mutate(neurodivergent =if_else(k2q35a ==1| k2q38a ==1| k2q36a ==1| k2q60a ==1| k2q37a ==1| k2q30a ==1| downsyn ==1| k2q31a ==1| k2q61a ==1| k2q42a ==1, 1, 0)) %>%# Filter data to include only individuals aged 6 or above and identified as neurodivergentfilter(sc_age_years >=6, neurodivergent ==1)# Resulting filtered data to be used for further analysiswrite.csv(neurodivergent_data_23, "neurodivergent_data_2023.csv")## Manually removed index column in both 2016-2022 and 2023 files ## Manually add empty `family` & `menbevcov` columns to 2023 file# Read in update 2023 data fileneurodivergent_data_23 <-read.csv("~/Downloads/nsch_2023e_topical_SAS/neurodivergent_data_2023.csv")# Read in prior years of neurodivergent dataneurodivergent_data_16to22 <-read.csv("~/Downloads/nsch_2023e_topical_SAS/neurodivergent_data_2016to2022.csv")# Merge two data frames with identical columnsneurodivergent_data_16to23 <-rbind(neurodivergent_data_16to22, neurodivergent_data_23)# Write the merged data frame into a CSV filewrite.csv(neurodivergent_data_16to23, "neurodivergent_data_2016to2023.csv", row.names =FALSE)
Analytic Data Preparation
Child Demographic Variables
Recode survey year into a factor, recode child sex as a categorical variable, recode race and ethnicity into descriptive categories, and create a combined race-ethnicity variable.
Expand Code
# Recode survey year into a factor variableneurodivergent_data <- neurodivergent_data %>%mutate(survey_year =factor(year))# Recode sex and race/ethnicity categoriesneurodivergent_data <- neurodivergent_data %>%mutate(# Recode sex into categorical factorssc_sex_cat =factor(sc_sex, levels =c(1, 2), labels =c("Male", "Female")),# Recode race categories into descriptive labelsrace =case_when( sc_race_r ==1~"White", sc_race_r ==2~"Black", sc_race_r ==3~"American Indian or Alaska Native", sc_race_r %in%4:5~"Asian, Native Hawaiian, or Pacific Islander", sc_race_r %in%6:7~"Other or Mixed Race",TRUE~NA_character_ ),# Recode ethnicityethnicity =case_when( sc_hispanic_r ==1~"Hispanic/Latino", sc_hispanic_r ==2~"Non-Hispanic/Latino",TRUE~NA_character_ ),# Create combined race and ethnicity categorysc_race_cat =case_when( race =="White"& ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic White", race =="Black"& ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic Black or African American", race =="American Indian or Alaska Native"& ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic American Indian or Alaska Native", race =="Asian, Native Hawaiian, or Pacific Islander"& ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic Asian, Native Hawaiian, or Pacific Islander", race =="Other or Mixed Race"& ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic Other or Mixed Race", ethnicity =="Hispanic/Latino"~"Hispanic/Latino of Any Race",TRUE~NA_character_ ),sc_race_cat =factor(sc_race_cat, levels =c("Non-Hispanic White", "Non-Hispanic Black or African American", "Non-Hispanic American Indian or Alaska Native", "Non-Hispanic Asian, Native Hawaiian, or Pacific Islander", "Non-Hispanic Other or Mixed Race", "Hispanic/Latino of Any Race" )) )
Healthcare-Related Variables
Create healthcare-related indices including Healthcare Access, Barriers to Access, Insurance Type, and Insurance Adequacy. Recode variables for consistent formatting, compute indices, and categorize healthcare access levels based on defined thresholds.
Expand Code
# Create a Healthcare Access Index and categorize itneurodivergent_data <- neurodivergent_data %>%mutate(healthcare_access_index =case_when( k4q27 ==2& k4q22_r ==1~2, # High access k4q27 ==2& k4q22_r ==3~1, # Medium access k4q27 ==1| k4q22_r ==2~0, # Low accessTRUE~NA_real_# Assign NA for missing or undefined values ),HAI_category =factor(healthcare_access_index, levels =c(2, 1, 0), labels =c("High Access", "Medium Access", "Low Access")) )# Recode 'treatneed' to standardize across survey yearsneurodivergent_data <- neurodivergent_data %>%mutate(treatneed_score =case_when( survey_year %in%2018:2023& treatneed ==1~0, # Not difficult (2018-2023) survey_year %in%2018:2023& treatneed ==2~1, # Somewhat difficult (2018-2023) survey_year %in%2018:2023& treatneed ==3~2, # Very difficult (2018-2023) survey_year %in%2018:2023& treatneed ==4~3, # Not possible (2018-2023) survey_year %in%2016:2017& treatneed ==1~0, # Not a problem (2016-2017) survey_year %in%2016:2017& treatneed ==2~1, # Small problem (2016-2017) survey_year %in%2016:2017& treatneed ==3~2, # Big problem (2016-2017)TRUE~NA_real_# Assign NA for missing or undefined values ) )# Create the Barriers to Access Indexneurodivergent_data <- neurodivergent_data %>%mutate(barriers_score =rowSums(select(., appointment, available, issuecost, notelig, notopen, transportcc, treatneed_score) ==1, na.rm =TRUE),barriers_cat =cut(barriers_score, breaks =c(-1, 0, 2, 4, 7), labels =c("No Barriers", "Few Barriers", "Moderate Barriers", "Many Barriers")) )# Include insurance type as a categorical variableneurodivergent_data <- neurodivergent_data %>%mutate(insurance_type =case_when( instype ==1~"Public Only", instype ==2~"Private Only", instype ==3~"Private + Public", instype ==5~"Uninsured", instype ==4~NA_character_# Excluded 'insurance type unspecified', n = 85 ),insurance_type =factor(insurance_type, levels =c("Private Only", "Private + Public", "Public Only", "Uninsured")) )# Create an Insurance Adequacy Indexneurodivergent_data <- neurodivergent_data %>%mutate(# Assign scores to the categorical responsesk3q20_score =recode(k3q20, `1`=4, `2`=3, `3`=2, `4`=1), k3q22_score =recode(k3q22, `1`=4, `2`=3, `3`=2, `4`=1),menbevcov_score =recode(menbevcov, `1`=4, `2`=3, `3`=2, `4`=1, `5`=NA_real_), # Excluded `menbevcov` = 5 'child does not need/use M/BH services', n = 15,310# Sum scores to create an Insurance Adequacy Indexinsurance_adequacy_index =rowSums(across(c(k3q20_score, k3q22_score, menbevcov_score)), na.rm =TRUE),# Categorize the Insurance Adequacy IndexIAI_category =case_when( insurance_adequacy_index >=10& insurance_adequacy_index <=12~"High Adequacy", insurance_adequacy_index >=7& insurance_adequacy_index <10~"Moderate Adequacy", insurance_adequacy_index >=3& insurance_adequacy_index <7~"Low Adequacy",TRUE~NA_character_# Assign NA for missing or undefined values ),IAI_category =factor(IAI_category, levels =c("High Adequacy", "Moderate Adequacy", "Low Adequacy")) )
1
k4q27 = Needed Health Care Not Received; k4q22_r = Mental Health Professional Treatment
2
treatneed = Mental Health Professional Treatment - Access Difficulty; appointment = Needed Health Care Not Received Due to - Getting Appointment; available = Needed Health Care Not Received Due to - Not Available; issuecost = Needed Health Care Not Received Due to - Cost; notelig = Needed Health Care Not Received Due to - Not Eligible; notopen = Needed Health Care Not Received Due to - Office Not Open; transportcc = Needed Health Care Not Received Due to - Getting Transportation
3
k3q20 = Benefits Cover Services; k3q22 = Allow to See Provider; menbevcov = Cover Mental Behavioral Needs
Family & Caregiver Context Variables
Recode caregiver education, family income, and caregiver composition variables to categorical values. Create a composite score for caregiver mental health, then categorize it into meaningful health levels.
Expand Code
# Recode caregiver education and family income categoriesneurodivergent_data <- neurodivergent_data %>%mutate(# Recode caregiver education into categoriescaregiver_education =case_when( higrade_tvis ==1~"Less Than High School", higrade_tvis ==2~"High School/Vocational", higrade_tvis ==3~"Some College/Associate Degree", higrade_tvis ==4~"College Degree/Higher",TRUE~NA_character_# Assign NA for missing or undefined values ),caregiver_education =factor(caregiver_education, levels =c("Less Than High School", "High School/Vocational", "Some College/Associate Degree", "College Degree/Higher")),# Recode family income into federal poverty level categoriesfpl_category =case_when( fpl %in%c(50:99) | fpl_mean <100~"<100%", fpl %in%c(100:199) | fpl_mean <200~"100%-199%", fpl %in%c(200:299) | fpl_mean <300~"200%-299%", fpl %in%c(300:399) | fpl_mean <400~"300%-399%", fpl %in%c(400:999) | fpl_mean <999~"≥400%",TRUE~NA_character_# Assign NA for missing or undefined values ),fpl_category =factor(fpl_category, levels =c("<100%", "100%-199%", "200%-299%", "300%-399%", "≥400%")) )# Simplify caregiver composition based on family structureneurodivergent_data <- neurodivergent_data %>%mutate(caregiver_composition =case_when( family %in%c(1, 2, 3, 4) | family_r %in%c(1, 2, 3, 4) ~"Dual-caregiver household", family %in%c(5, 6, 7, 8, 9) | family_r %in%c(5, 6, 7, 8) ~"Single-caregiver household",TRUE~NA_character_# Assign NA for missing or undefined values ),caregiver_composition =factor(caregiver_composition, levels =c("Dual-caregiver household", "Single-caregiver household")) )# Create a composite variable for caregiver mental healthneurodivergent_data <- neurodivergent_data %>%mutate(caregiver_mental_health_score =rowMeans(select(., a1_menthealth, a2_menthealth), na.rm =TRUE),caregiver_mental_health =cut( caregiver_mental_health_score,breaks =c(1, 1.5, 2.5, 3.5, 4.5, 5),labels =c("Excellent", "Very Good", "Good", "Fair", "Poor"),include.lowest =TRUE ) )
1
Scores of caregiver mental health are averaged to account for single respondents (i.e., single-caregiver households).
Neighborhood Context Variables
Recode neighborhood binary items to 0/1 format and standardize ordinal neighborhood variables by subtracting their mean and dividing by the standard deviation.
k10q11 = Neighborhood - Sidewalks or Walking Paths; k10q12 = Neighborhood - Park or Playground; k10q13 = Neighborhood - Recreation Center; k10q14 = Neighborhood - Library or Bookmobile; k10q20 = Neighborhood - Litter or Garbage; k10q22 = Neighborhood - Poorly Kept or Rundown Housing; k10q23 = Neighborhood - Vandalism; k9q96 = Other Adult Child Can Rely On For Advice
2
k10q40_r = Child is Safe In Neighborhood; k10q41_r = Child Is Safe at School; k10q30 = People In Neighborhood Help Each Other Out; k10q31 = Watch Out for Other’s Children
Factor Analysis - Neighborhood Context Scores
Conduct factor analysis on neighborhood context variables to derive factor scores for physical resources, environmental threats, and collective efficacy, and add these scores to the primary data frame. Remove intermediate data frames from the global environment.
The factor structure was adapted from factor conceptualization detailed in Fan & Chen (2012), “Family functioning as a mediator between neighborhood conditions and children’s health: Evidence from a national survey in the United States.” Principal Axis Factoring (PAF) was used as the factoring method, and factor scores were calculated using the regression method, thus consistent with the approach in the referenced study.
Expand Code
# Select the items for each neighborhood factorphysical_resources_items <- neurodivergent_data %>%select(k10q11, k10q12, k10q13, k10q14)environmental_threats_items <- neurodivergent_data %>%select(k10q20, k10q22, k10q23, k10q40_r, k10q41_r)collective_efficacy_items <- neurodivergent_data %>%select(k10q40_r, k10q41_r, k10q30, k10q31, k9q96)# Conduct factor analysis for each factor using Principal Axis Factoring (PAF)fa_physical_resources <-fa(physical_resources_items, nfactors =1, fm ="pa")fa_environmental_threats <-fa(environmental_threats_items, nfactors =1, fm ="pa")fa_collective_efficacy <-fa(collective_efficacy_items, nfactors =1, fm ="pa")# Calculate factor scores for physical resources using regression methodphysical_resources_scores <-factor.scores(physical_resources_items, fa_physical_resources, method ="regression")$scores[, 1]# Calculate factor scores for environmental threats using regression methodenvironmental_threats_scores <-factor.scores(environmental_threats_items, fa_environmental_threats, method ="regression")$scores[, 1]# Calculate factor scores for collective efficacy using regression methodcollective_efficacy_scores <-factor.scores(collective_efficacy_items, fa_collective_efficacy, method ="regression")$scores[, 1]# Add factor scores back to the datasetneurodivergent_data <- neurodivergent_data %>%mutate(physical_resources = physical_resources_scores,environmental_threats = environmental_threats_scores,collective_efficacy = collective_efficacy_scores )# Remove intermediate data frames from the global environmentrm(physical_resources_items, environmental_threats_items, collective_efficacy_items, fa_physical_resources, fa_environmental_threats, fa_collective_efficacy, physical_resources_scores, environmental_threats_scores, collective_efficacy_scores)
Child Health & Development Variables
Standardize Down Syndrome & Tourette Syndrome severity variables based on survey year, calculate total severity scores across diagnoses, categorize average severity scores into severity levels, recode child general health, and convert the special healthcare needs indicator into binary format with labeled categories.
Expand Code
# Standardize Down Syndrome and Tourette Syndrome severity variables based on survey yearneurodivergent_data <- neurodivergent_data %>%mutate(downsyn_desc_recode =case_when( survey_year %in%c(2016, 2017) & downsyn_desc ==1~1, survey_year %in%c(2016, 2017) & downsyn_desc ==4~3, survey_year >=2019~ downsyn_desc,TRUE~NA_real_# Handle missing or undefined years as NA ),k2q38c_recode =case_when( k2q38c ==1~1, k2q38c ==4~3 ) )# List of severity variables, including recoded Down & Tourette Syndromes' descriptionseverity_variables <-c("k2q31c", "k2q35c", "k2q38c_recode", "k2q30c", "k2q36c", "k2q60c", "k2q37c", "downsyn_desc_recode", "cerpals_desc", "k2q42c")# Calculate total and average severity scoresneurodivergent_data <- neurodivergent_data %>%rowwise() %>%mutate(total_severity_score =ifelse(all(is.na(c_across(all_of(severity_variables)))), NA_real_, sum(c_across(all_of(severity_variables)), na.rm =TRUE)) ) %>%ungroup()# Recode child general health variableneurodivergent_data <- neurodivergent_data %>%mutate(child_gen_health =factor(case_when( k2q01 ==1~"Excellent", k2q01 ==2~"Very Good", k2q01 ==3~"Good", k2q01 ==4~"Fair", k2q01 ==5~"Poor" ),levels =c("Excellent", "Very Good", "Good", "Fair", "Poor")) )# Convert Special Healthcare Needs (SHCN) indicator to binary and label accordinglyneurodivergent_data <- neurodivergent_data %>%mutate(sc_cshcn =ifelse(sc_cshcn ==1, 1, 0),special_healthcare_needs =factor(case_when( sc_cshcn ==0~"Non-SHCN", sc_cshcn ==1~"SHCN" ),levels =c("Non-SHCN", "SHCN") ) )
Recode ACEs, psychopathology, and psychosocial factors. Compute the Childhood Flourishing Index and the Family Resilience and Connection Index, based on Bethell et al. (2019), “Family Resilience And Connection Promote Flourishing Among US Children, Even Amid Adversity.” Create categorical versions of these indices for further analysis.
The Child Flourishing Index (CFI) in was calculated using three items from the National Survey of Children’s Health (NSCH). Each item asked parents how well each of the following statements described their child: (1) “shows interest and curiosity in learning new things,” (2) “works to finish tasks he or she starts,” and (3) “stays calm and in control when faced with a challenge.” Each item was assigned one point for a response of “definitely true.” Children who scored 3 were classified as flourishing.
The Family Resilience and Connection Index (FRCI) was created using six items from the NSCH. The index included four items on family resilience, which asked how often families: (1) “talk together about what to do,” (2) “work together to solve our problems,” (3) “know we have strengths to draw on,” and (4) “stay hopeful even in difficult times.” Parents could earn one point for each of these items if they responded “all of the time.” Two additional items assessed family connection: (1) the ability of the parent to “share ideas or talk about things that really matter” with their child and (2) how well parents were “handling the day-to-day demands of raising children.” One point was assigned for each response of “very well.” The FRCI score ranged from 0 to 6.
Expand Code
# Recode Adverse Childhood Experiences (ACE) variables & psychopathology outcomesneurodivergent_data <- neurodivergent_data %>%mutate(# Recode 'ace1' to a dichotomous variableace1_recode =ifelse(ace1 %in%c(2, 3, 4), 1, ifelse(ace1 ==1, 0, NA)) ) %>%mutate(# Total ACE count excluding missing valuesace_total =rowSums(select(., ace1_recode, ace3:ace12) ==1, na.rm =TRUE),# Categorize total ACEs for analysisace_counts =cut(ace_total, breaks =c(-1, 0, 1, 3, Inf), labels =c('0 ACEs', '1 ACE', '2-3 ACEs', '4+ ACEs'), right =TRUE),ace_counts =factor(ace_counts, levels =c("0 ACEs", "1 ACE", "2-3 ACEs", "4+ ACEs")),# Recode psychopathology outcomesAnxiety =factor(k2q33b, levels =c(2, 1), labels =c("No", "Yes")),Depression =factor(k2q32b, levels =c(2, 1), labels =c("No", "Yes")),Behavioral_Problems =factor(k2q34b, levels =c(2, 1), labels =c("No", "Yes")) )# Calculate and recode Childhood Flourishing Index (CFI)neurodivergent_data <- neurodivergent_data %>%mutate(# Calculate the Childhood Flourishing Index (CFI)CFI =rowSums(select(., k6q71_r, k7q85_r, k7q84_r) ==1, na.rm =TRUE),# Recode CFI into a dichotomous variableCFI_dich =case_when( CFI ==3~"Flourishing", CFI %in%0:2~"Not Flourishing",TRUE~NA_character_ ),CFI_dich =factor(CFI_dich, levels =c("Not Flourishing", "Flourishing")) )# Calculate and recode Family Resilience and Connection Index (FRCI)neurodivergent_data <- neurodivergent_data %>%mutate(# Recode FRCI items: 1 point for "All of the time" (response = 1), handle missing valuestalkabout_score =ifelse(!is.na(talkabout) & talkabout ==1, 1, 0),wktosolve_score =ifelse(!is.na(wktosolve) & wktosolve ==1, 1, 0),strengths_score =ifelse(!is.na(strengths) & strengths ==1, 1, 0),hopeful_score =ifelse(!is.na(hopeful) & hopeful ==1, 1, 0),# Recode additional FRCI items: 1 point for "Very well" (response = 1), handle missing valuesk8q21_score =ifelse(!is.na(k8q21) & k8q21 ==1, 1, 0),k8q30_score =ifelse(!is.na(k8q30) & k8q30 ==1, 1, 0),# Calculate total FRCI score by summing all items, ignoring missing valuesFRCI =rowSums(across(c(talkabout_score, wktosolve_score, strengths_score, hopeful_score, k8q21_score, k8q30_score)), na.rm =TRUE),# Create categorical version of FRCIFRCI_cat =factor(case_when( FRCI <=1~"0-1", FRCI <=3~"2-3", FRCI <=6~"4-6",TRUE~NA_character_ ),levels =c("0-1", "2-3", "4-6") ) )
1
ace1 is originally coded by NSCH as a frequency measure, but it’s needed as indicative variable of presence like the other ace# variables.
2
ace1 = Hard to Cover Basics Like Food or Housing; ace3 = Child Experienced - Parent or Guardian Divorced; ace4 = Child Experienced - Parent or Guardian Died; ace5 = Child Experienced - Parent or Guardian Time in Jail; ace6 = Child Experienced - Adults Slap, Hit, Kick, Punch Others; ace7 = Child Experienced - Victim of Violence; ace8 = Child Experienced - Lived with Mentally Ill; ace9 = Child Experienced - Lived with Person with Alcohol/Drug Problem; ace10 = Child Experienced - Treated Unfairly Because of Race; ace11 = Child Experienced - Treated Unfairly Because of Health Condition; ace12 = Child Experienced - Treated Unfairly Because of their Sexual Orientation or Gender Identity
k6q71_r = Show Interest and Curiosity; k7q85_r = Stays Calm and In Control When Challenged; k7q84_r = Works to Finish Tasks Started
5
talkabout = Facing Problems - How Often Talk Together; wktosolve = Facing Problems - How Often Work Together; strengths = Facing Problems - How Often Draw on Strengths; k8q21 = Share Ideas or Talk About Things That Matter; k8q30 = How Well Handling Demands of Raising Children
Frequencies & Descriptive Statistics
Sample Trends of Neurodevelopmental Diagnoses
Calculate yearly prevalence rates for neurodevelopmental disorder (NDD) diagnoses without survey weights. Group data by survey year, calculate prevalence for each diagnosis, and transform it into a long format. Create a line plot to visualize yearly trends in the prevalence of each NDD diagnosis, using color to differentiate between the diagnoses.
Expand Code
# Generate list of NDD diagnoses among subject samplevariables <-c("k2q31a", "k2q35a", "k2q38a", "k2q30a", "k2q36a", "k2q60a", "k2q37a", "downsyn", "k2q61a", "k2q42a")variable_names <-c("ADHD", "Autism/ASD", "Tourette’s Syndrome", "Learning Disability","Development Delay", "Intellectual Disability", "Speech/Other Language Disorder", "Down Syndrome", "Cerebral Palsy", "Epilepsy")# Calculate yearly prevalence rates without survey weights for specific samplesample_prevalence <- neurodivergent_data %>%group_by(survey_year) %>%summarise(across(all_of(variables), ~mean(. ==1, na.rm =TRUE) *100, .names ="Prevalence_{col}")) %>%pivot_longer(cols =starts_with("Prevalence_"), names_to ="Diagnosis", values_to ="Prevalence" ) %>%mutate(Diagnosis =gsub("Prevalence_", "", Diagnosis),Diagnosis =factor(Diagnosis, levels = variables, labels = variable_names) )# Plot the prevalence rates by year for each NDD diagnosis without survey weightsggplot(sample_prevalence, aes(x = survey_year, y = Prevalence, group = Diagnosis, color = Diagnosis)) +geom_line(size =1.1) +geom_point(size =2) +labs(title ="Yearly Prevalence of Neurodevelopmental Disorders (Sample Estimates)",x ="Survey Year",y ="Prevalence (%)",color ="Diagnosis" ) +theme_minimal() +theme(axis.text.x =element_text(angle =45, hjust =1))
1
The following diagnoses’ variables correspond to whether a child has EVER received such a diagnosis: k2q31a = ADD/ADHD; k2q35a = Autism ASD; k2q38a = Tourette Syndrome; k2q30a = Learning Disability; k2q36a = Developmental Delay; k2q60a = Intellectual Disability; k2q37a = Speech Disorder; downsyn = Down Syndrome; cerpals = Cerebral Palsy; k2q42a = Epilepsy
Expand Code
# Remove intermediate data frames from the global environmentrm(design, sample_prevalence, yearly_prevalence, yearly_prevalence_df)
Sample Demographic Distributions
Calculate and visualize the distribution of child age, educational attainment, sex, socioeconomic status, and race among neurodivergent participants. Create a raincloud plot for age distribution alongside a summary statistics table, and generate bar plots for key categorical variables, ensuring all visuals are clear and informative.
Expand Code
# Calculate descriptives of age for neurodivergent participantsneurodiv_age_summary <- neurodivergent_data %>%summarise(Mean =mean(sc_age_years, na.rm =TRUE),n =n(),SE =sd(sc_age_years, na.rm =TRUE) /sqrt(n()),SD =sd(sc_age_years, na.rm =TRUE),Minimum =min(sc_age_years, na.rm =TRUE),Q1 =quantile(sc_age_years, 0.25, na.rm =TRUE),Median =median(sc_age_years, na.rm =TRUE),Q3 =quantile(sc_age_years, 0.75, na.rm =TRUE),Maximum =max(sc_age_years, na.rm =TRUE),IQR =quantile(sc_age_years, 0.75, na.rm =TRUE) -quantile(sc_age_years, 0.25, na.rm =TRUE),.groups ='drop' )# Create raincloud plot displaying child age descriptivesraincloud_plot <- neurodivergent_data %>%ggplot(aes(x =1, y = sc_age_years)) + PupillometryR::geom_flat_violin(adjust =1.5, trim =FALSE, alpha =0.5, color =NA, fill ="deeppink", position =position_nudge(x =0.2) ) +geom_boxplot(width =0.12, outlier.color ="black", colour ="black", fill =NA, alpha =0.5, position =position_nudge(x =0.2) ) +stat_dots(dotsize =0.1, side ="left", justification =1, binwidth =1, alpha =0.5, colour ="deeppink", fill ="deeppink" ) +geom_point(data = neurodiv_age_summary, aes(x =1, y = Mean), shape =18, color ="deeppink", size =3, position =position_nudge(x =0.2)) +geom_errorbar(data = neurodiv_age_summary, aes(x =1, y = Mean, ymin = Mean - SE, ymax = Mean + SE), width =0.1, color ="deeppink", position =position_nudge(x =0.2)) +labs(title ="Age Distribution of Neurodivergent Participants", y ="Age (Years)", x ="") +scale_y_continuous(breaks =seq(5, 18, by =1), limits =c(5, 18)) +coord_flip() +theme_minimal() +theme(axis.text.y =element_blank(), legend.position ="none")# Create summary statistics table as a data frame and transpose itsummary_df <- neurodiv_age_summary %>%select(Mean, SD, SE, Median, Minimum, Maximum, IQR, Q1, Q3) %>%t() %>%as.data.frame() %>%mutate_all(~round(., 2))# Update row names to meaningful labelscolnames(summary_df) <-"Value"summary_df$Statistic <-rownames(summary_df)summary_df <- summary_df %>%select(Statistic, Value)# Convert the transposed summary statistics data frame to a table using gridExtrasummary_table <- gridExtra::tableGrob(summary_df, rows =NULL)# Combine the plot and summary table using gridExtragridExtra::grid.arrange(raincloud_plot, summary_table, nrow =1)
# Bar plot for sex counts distributionggplot(neurodivergent_data, aes(x = sc_sex_cat)) +geom_bar(aes(fill = sc_sex_cat), show.legend =FALSE) +geom_text(stat ='count', aes(label = ..count..), position =position_stack(vjust =0.5), color ="white") +scale_fill_manual(values =c("Male"="lightblue", "Female"="pink")) +labs(title ="Sex Distribution of Neurodivergent Participants",x ="Sex",y ="Count" ) +theme_minimal()
Expand Code
# Bar plot for SES counts categoriesggplot(neurodivergent_data, aes(x = fpl_category)) +geom_bar(aes(fill = fpl_category), show.legend =FALSE) +geom_text(stat ='count', aes(label = ..count..), position =position_stack(vjust =0.5), color ="white") +labs(title ="Distribution of Socioeconomic Status (FPL)",x ="Federal Poverty Level Groups",y ="Count" ) +theme_minimal()
Expand Code
# Bar plot for race counts categoriesggplot(neurodivergent_data, aes(x = sc_race_cat)) +geom_bar(aes(fill = sc_race_cat), show.legend =FALSE) +geom_text(stat ='count', aes(label = ..count..), vjust =-0.3, size =3.5) +labs(title ="Race Distribution of Neurodivergent Participants",x ="Race",y ="Count" ) +theme_minimal() +theme(axis.text.x =element_text(angle =45, hjust =1))
Expand Code
# Remove intermediate data frames from the global environmentrm(neurodiv_age_summary, raincloud_plot, summary_df, summary_table)
Sample Trends of Child Flourishing
Create visualizations for the Childhood Flourishing Index (CFI) and dichotomous CFI variable (Flourishing vs Not Flourishing). Generate bar plots to show the distribution of these variables and line plots to display yearly prevalence rates across survey years.
Expand Code
# Bar plot for Childhood Flourishing Index (CFI) Score distributionggplot(neurodivergent_data, aes(x =factor(CFI), fill =factor(CFI))) +geom_bar(show.legend =FALSE) +geom_text(stat ='count', aes(label = ..count..), position =position_stack(vjust =0.5), color ="white") +labs(x ="Childhood Flourishing Index (CFI) Score",y ="Number of Children",title ="Distribution of Children by CFI Score" ) +theme_minimal() +theme(axis.text.x =element_text(angle =0))
Expand Code
# Calculate prevalence of CFI scores grouped by yearflourishing_prevalence <- neurodivergent_data %>%group_by(survey_year, CFI) %>%summarise(Count =n(),.groups ="drop" ) %>%group_by(survey_year) %>%mutate(Total =sum(Count),Prevalence = (Count / Total) *100 )# Line plot for yearly prevalence of CFI scoresggplot(flourishing_prevalence, aes(x = survey_year, y = Prevalence, color =as.factor(CFI), group = CFI)) +geom_line(size =1.2) +geom_point(size =3) +labs(title ="Yearly Prevalence of Child Flourishing Index Scores (Sample Estimate)",x ="Survey Year",y ="Prevalence (%)",color ="Child Flourishing Index Score" ) +theme_minimal()
Expand Code
# Bar plot for CFI dichotomous variable distributionggplot(neurodivergent_data, aes(x = CFI_dich, fill = CFI_dich)) +geom_bar(show.legend =FALSE) +geom_text(stat ='count', aes(label = ..count..), position =position_stack(vjust =0.5), color ="white") +labs(title ="Distribution of CFI Dichotomous Variable",x ="Child Flourishing Index (Dichotomous)",y ="Count" ) +theme_minimal() +theme(axis.text.x =element_text(angle =0))
Expand Code
# Calculate prevalence of flourishing status (CFI_dich) grouped by yearflourishing_prevalence_status <- neurodivergent_data %>%group_by(survey_year, CFI_dich) %>%summarise(Count =n(),.groups ="drop" ) %>%group_by(survey_year) %>%mutate(Total =sum(Count),Prevalence = (Count / Total) *100 )# Line plot for yearly prevalence of flourishing statusggplot(flourishing_prevalence_status, aes(x = survey_year, y = Prevalence, color = CFI_dich, group = CFI_dich)) +geom_line(size =1.2) +geom_point(size =3) +labs(title ="Yearly Prevalence of Child Flourishing Status (Sample Estimate)",x ="Survey Year",y ="Prevalence (%)",color ="Child Flourishing Status" ) +theme_minimal()
Expand Code
# Remove intermediate data frames from the global environmentrm(flourishing_prevalence, flourishing_prevalence_status)
Sample Trends of Family Resilience & Connection
Create visualizations for the Family Resilience and Connection Index (FRCI) scores and categories. Generate bar plots to show the distribution of FRCI scores and categories, and line plots to display yearly prevalence rates across survey years for both continuous FRCI scores and categorical FRCI status.
Expand Code
# Bar plot for Family Resilience and Connection Index (FRCI) Score distributionggplot(neurodivergent_data, aes(x =factor(FRCI), fill =factor(FRCI))) +geom_bar(show.legend =FALSE) +geom_text(stat ='count', aes(label = ..count..), position =position_stack(vjust =0.5), color ="white") +labs(x ="Family Resilience and Connection Index (FRCI) Score",y ="Number of Children",title ="Distribution of Children by FRCI Score" ) +theme_minimal() +theme(axis.text.x =element_text(angle =0))
Expand Code
# Calculate prevalence of FRCI scores grouped by yearfrci_prevalence_scores <- neurodivergent_data %>%group_by(survey_year, FRCI) %>%summarise(Count =n(),.groups ="drop" ) %>%group_by(survey_year) %>%mutate(Total =sum(Count),Prevalence = (Count / Total) *100 )# Line plot for yearly prevalence of FRCI scoresggplot(frci_prevalence_scores, aes(x = survey_year, y = Prevalence, color =as.factor(FRCI), group = FRCI)) +geom_line(size =1.2) +geom_point(size =3) +labs(title ="Yearly Prevalence of Family Resilience & Connection Index (FRCI) Scores (Sample Estimate)",x ="Survey Year",y ="Prevalence (%)",color ="FRCI Score" ) +theme_minimal()
Expand Code
# Bar plot for FRCI categorical variable distributionggplot(neurodivergent_data, aes(x = FRCI_cat, fill = FRCI_cat)) +geom_bar(show.legend =FALSE) +geom_text(stat ='count', aes(label = ..count..), position =position_stack(vjust =0.5), color ="white") +labs(title ="Distribution of FRCI Categorical Variable",x ="Family Resilience and Connection Index (Categorical)",y ="Count" ) +theme_minimal() +theme(axis.text.x =element_text(angle =0))
Expand Code
# Calculate prevalence of FRCI categories grouped by yearfrci_prevalence_categories <- neurodivergent_data %>%group_by(survey_year, FRCI_cat) %>%summarise(Count =n(),.groups ="drop" ) %>%group_by(survey_year) %>%mutate(Total =sum(Count),Prevalence = (Count / Total) *100 )# Line plot for yearly prevalence of FRCI categoriesggplot(frci_prevalence_categories, aes(x = survey_year, y = Prevalence, color = FRCI_cat, group = FRCI_cat)) +geom_line(size =1.2) +geom_point(size =3) +labs(title ="Yearly Prevalence of Family Resilience & Connection Index Categories (Sample Estimate)",x ="Survey Year",y ="Prevalence (%)",color ="FRCI Category" ) +theme_minimal()
Expand Code
# Remove intermediate data frames from the global environmentrm(frci_prevalence_scores, frci_prevalence_categories)
Sample Trends of Psychopathology
Create a stacked bar plot to show the distribution of psychopathology outcomes across participants. Calculate the yearly prevalence of anxiety, depression, and behavioral problems, and visualize these trends using a line plot.
Expand Code
# Create a long-format data frame for plotting psychopathology outcomespsychopathology_data <- neurodivergent_data %>%select(Anxiety, Depression, Behavioral_Problems) %>%pivot_longer(cols =everything(), names_to ="Psychopathology", values_to ="Response") %>%filter(!is.na(Response)) # Exclude missing responses# Create a stacked bar plot for psychopathology outcome distributionggplot(psychopathology_data, aes(x = Psychopathology, fill = Response)) +geom_bar(position ="stack") +geom_text(stat ='count', aes(label = ..count..), position =position_stack(vjust =0.5), color ="white") +scale_fill_manual(values =c("Yes"="blue", "No"="red")) +labs(title ="Distribution of Psychopathology Outcomes",x ="Psychopathology Outcome",y ="Count",fill ="Response" ) +theme_minimal()
Expand Code
# List of psychopathology outcomes for analysispsychopathology_variables <-c("Anxiety", "Depression", "Behavioral_Problems")# Calculate prevalence of each psychopathology outcome by survey yearpsychopathology_prevalence <-lapply(psychopathology_variables, function(var) { neurodivergent_data %>%group_by(survey_year, !!sym(var)) %>%summarise(Count =n(), # Count of each level.groups ="drop" ) %>%group_by(survey_year) %>%mutate(Total =sum(Count), # Total count per yearPrevalence = (Count / Total) *100, # Calculate prevalence as a percentageOutcome = var # Add the outcome name ) %>%filter(!!sym(var) =="Yes") # Only consider "Yes" (those with the outcome)}) %>%bind_rows() # Combine all results into one data frame# Visualize the prevalence of psychopathology outcomes by year using a line plotggplot(psychopathology_prevalence, aes(x = survey_year, y = Prevalence, color = Outcome, group = Outcome)) +geom_line(size =1.2) +# Add linesgeom_point(size =3) +# Add pointslabs(title ="Yearly Prevalence of Psychopathology Outcomes (Sample Estimate)",x ="Survey Year",y ="Prevalence (%)",color ="Outcome" ) +theme_minimal() +theme(axis.text.x =element_text(angle =0) )
Expand Code
# Remove intermediate data frames from the global environmentrm(psychopathology_data, psychopathology_prevalence)
Distribution & Cross-Tabulations of ACEs
Create a bar plot to visualize the distribution of children by ACE counts. Conduct cross-tabulations between ACE counts and psychopathology/flourishing variables, as well as between psychopathology outcomes and the flourishing variable, displaying the results in table format.
Expand Code
# Bar plot of ACE counts categoriesggplot(neurodivergent_data, aes(x =factor(ace_counts), fill =factor(ace_counts))) +geom_bar(show.legend =FALSE) +geom_text(stat ='count', aes(label = ..count..), position =position_stack(vjust =0.5), color ="white") +labs(x ="ACE Counts", y ="Number of Children", fill ="ACE Counts",title ="Distribution of Children by ACE Counts") +theme_minimal() +theme(axis.text.x =element_text(angle =0))
1
Cross-tabulations (or contingency tables) summarize the relationship between two or more categorical variables. In this case, I’m exploring how the mental health conditions (Anxiety, Depression, Behavioral Issues) and Child Flourishing Index (CFI) are distributed across different levels of ACEs.
Expand Code
# Cross-tabulation between psychopathology/flourishing variables and ACEs count variablecrosstable(neurodivergent_data, c(Anxiety, Depression, Behavioral_Problems, CFI_dich), by = ace_counts, showNA ="no", percent_digits =1, total ="both") %>%as_flextable()# Cross-tabulation between psychopathology outcomes and flourishing variablecrosstable(neurodivergent_data, c(Anxiety, Depression, Behavioral_Problems), by = CFI_dich, showNA ="no", percent_digits =1, total ="both") %>%as_flextable()
label
variable
ace_counts
Total
0 ACEs
1 ACE
2-3 ACEs
4+ ACEs
Anxiety
No
383 (27.0%)
375 (26.5%)
407 (28.7%)
252 (17.8%)
1417 (7.7%)
Yes
3389 (20.0%)
4366 (25.8%)
5180 (30.6%)
3992 (23.6%)
16927 (92.3%)
Total
3772 (20.6%)
4741 (25.8%)
5587 (30.5%)
4244 (23.1%)
18344 (100.0%)
Depression
No
310 (20.3%)
349 (22.8%)
486 (31.8%)
383 (25.1%)
1528 (16.6%)
Yes
1004 (13.1%)
1475 (19.3%)
2553 (33.4%)
2619 (34.2%)
7651 (83.4%)
Total
1314 (14.3%)
1824 (19.9%)
3039 (33.1%)
3002 (32.7%)
9179 (100.0%)
Behavioral_Problems
No
592 (21.7%)
692 (25.4%)
875 (32.1%)
569 (20.9%)
2728 (14.5%)
Yes
3052 (18.9%)
4234 (26.3%)
4939 (30.6%)
3890 (24.1%)
16115 (85.5%)
Total
3644 (19.3%)
4926 (26.1%)
5814 (30.9%)
4459 (23.7%)
18843 (100.0%)
CFI_dich
Not Flourishing
13914 (27.4%)
14990 (29.6%)
13868 (27.3%)
7936 (15.7%)
50708 (93.0%)
Flourishing
1530 (40.1%)
1239 (32.4%)
808 (21.2%)
243 (6.4%)
3820 (7.0%)
Total
15444 (28.3%)
16229 (29.8%)
14676 (26.9%)
8179 (15.0%)
54528 (100.0%)
label
variable
CFI_dich
Total
Not Flourishing
Flourishing
Anxiety
No
1332 (94.0%)
85 (6.0%)
1417 (7.7%)
Yes
16579 (97.9%)
348 (2.1%)
16927 (92.3%)
Total
17911 (97.6%)
433 (2.4%)
18344 (100.0%)
Depression
No
1477 (96.7%)
51 (3.3%)
1528 (16.6%)
Yes
7530 (98.4%)
121 (1.6%)
7651 (83.4%)
Total
9007 (98.1%)
172 (1.9%)
9179 (100.0%)
Behavioral_Problems
No
2633 (96.5%)
95 (3.5%)
2728 (14.5%)
Yes
15943 (98.9%)
172 (1.1%)
16115 (85.5%)
Total
18576 (98.6%)
267 (1.4%)
18843 (100.0%)
Inferential Statistics
Interaction Models
Fit three logistic regression models with survey weights to examine the moderating effect of Child Flourishing Index (CFI) on the association between Adverse Childhood Experiences (ACEs) and psychopathology outcomes (anxiety, depression, behavioral problems). Extract residuals, check heteroscedasticity using Breusch-Pagan tests, calculate Cook’s distance for influential observations, and obtain robust standard errors. Create summary tables with robust standard errors and False Discovery Rate (FDR) adjusted p-values for each model.
# Create survey design object with Final Weight for Child (FWC) as weightsdesign <-svydesign(ids =~1, weights =~fwc, data = neurodivergent_data) # will be used across all models# Logistic regression model with CFI_dich as a moderator for Anxiety, using robust standard errorsmodel_anx_mod <-svyglm(Anxiety ~ ace_counts * CFI_dich + sc_race_cat + sc_sex_cat + fpl_category + caregiver_education + sc_age_years + survey_year + caregiver_mental_health + child_gen_health + caregiver_composition + physical_resources + environmental_threats + collective_efficacy + special_healthcare_needs + total_severity_score + HAI_category + barriers_cat + insurance_type + IAI_category, design = design, family =quasibinomial())# Extract residuals from the logistic regression modelresiduals_anx <-residuals(model_anx_mod, type ="deviance")# Plot residuals to visually inspect themplot(residuals_anx, main ="Residuals Plot for Anxiety Model", ylab ="Residuals", xlab ="Index", pch =20, col ="blue")abline(h =0, col ="red", lty =2)
Expand Code
# Perform the Breusch-Pagan test for heteroscedasticitybptest(model_anx_mod)# Calculate Cook's distance for influential observations (leverage) in each modelcooks_anx <-cooks.distance(model_anx_mod)# Plot for Anxiety modelplot(cooks_anx, type ="p", main ="Cook's Distance for Anxiety Model",ylab ="Cook's Distance", xlab ="Observation Index")abline(h =4/length(cooks_anx), col ="red", lty =2) # Add a reference line for influential points
Expand Code
# Obtain robust standard errors for Anxiety modelrobust_anx <-coeftest(model_anx_mod, vcov =vcovHC(model_anx_mod, type ="HC0"))# Create a summary table for the Anxiety model with robust standard errorstbl_anx_mod <-tbl_regression(model_anx_mod, exponentiate =TRUE, label =list( ace_counts ~"ACEs", CFI_dich ~"Child Flourishing Index", fpl_category ~"Federal Poverty Level", caregiver_education ~"Caregiver Education", sc_age_years ~"Child Age", sc_race_cat ~"Child Race", sc_sex_cat ~"Child Sex", survey_year ~"Survey Year", caregiver_mental_health ~"Caregiver Mental Health", child_gen_health ~"Child General Health", caregiver_composition ~"Caregiver Composition", physical_resources ~"Neighborhood Physical Resources", environmental_threats ~"Neighborhood Environmental Threats", collective_efficacy ~"Neighborhood Collective Efficacy", special_healthcare_needs ~"Special Healthcare Needs", total_severity_score ~"NDD Severity (Total)", HAI_category ~"Healthcare Access Index", barriers_cat ~"Access Barriers Index", insurance_type ~"Insurance Type", IAI_category ~"Insurance Adequacy Index" )) %>%add_q() %>%# adjusts p-values for multiple testing (FDR)bold_p(t =0.05) %>%bold_p(t =0.05, q =TRUE) %>%# now bold q-values under the threshold of 0.05italicize_levels() %>% gtsummary::as_gt() %>% gt::tab_source_note(gt::md("*Note*. This logistic regression model examines the moderating effect of the Child Flourishing Index (CFI) on the association between Adverse Childhood Experiences (ACEs) and anxiety outcomes in neurodivergent children. Survey weights from the National Survey of Children's Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions and survey year, with FDR-adjusted p-values to account for multiple comparisons."))# Print table to review resultstbl_anx_mod
studentized Breusch-Pagan test
data: model_anx_mod
BP = 851.64, df = 52, p-value < 2.2e-16
Characteristic
OR1
95% CI1
p-value
q-value2
ACEs
0 ACEs
—
—
1 ACE
1.14
0.76, 1.71
0.5
0.8
2-3 ACEs
1.40
0.89, 2.18
0.14
0.3
4+ ACEs
1.55
0.95, 2.54
0.079
0.2
Child Flourishing Index
Not Flourishing
—
—
Flourishing
1.18
0.48, 2.89
0.7
0.9
Child Race
Non-Hispanic White
—
—
Non-Hispanic Black or African American
0.78
0.45, 1.36
0.4
0.6
Non-Hispanic American Indian or Alaska Native
0.57
0.23, 1.44
0.2
0.5
Non-Hispanic Asian, Native Hawaiian, or Pacific Islander
0.33
0.17, 0.65
0.001
0.012
Non-Hispanic Other or Mixed Race
0.94
0.61, 1.46
0.8
0.9
Hispanic/Latino of Any Race
0.61
0.43, 0.86
0.005
0.023
Child Sex
Male
—
—
Female
0.97
0.72, 1.31
0.9
0.9
Federal Poverty Level
<100%
—
—
100%-199%
0.58
0.34, 1.00
0.049
0.15
200%-299%
0.76
0.42, 1.36
0.4
0.6
300%-399%
0.51
0.28, 0.93
0.028
0.10
≥400%
0.75
0.42, 1.34
0.3
0.6
Caregiver Education
Less Than High School
—
—
High School/Vocational
3.11
1.42, 6.80
0.004
0.021
Some College/Associate Degree
2.99
1.44, 6.21
0.003
0.021
College Degree/Higher
4.25
2.00, 9.01
<0.001
0.003
Child Age
0.99
0.95, 1.03
0.7
0.9
Survey Year
2016
—
—
2017
1.14
0.61, 2.14
0.7
0.9
2018
0.92
0.56, 1.52
0.8
0.9
2019
0.65
0.40, 1.07
0.087
0.2
2020
1.21
0.68, 2.15
0.5
0.8
2021
1.90
1.22, 2.96
0.004
0.021
2022
1.44
0.89, 2.33
0.14
0.3
2023
1.24
0.76, 2.02
0.4
0.6
Caregiver Mental Health
Excellent
—
—
Very Good
1.25
0.90, 1.73
0.2
0.4
Good
1.78
1.21, 2.60
0.003
0.021
Fair
2.45
1.43, 4.20
0.001
0.012
Poor
5.86
1.74, 19.7
0.004
0.021
Child General Health
Excellent
—
—
Very Good
1.03
0.77, 1.38
0.9
0.9
Good
1.58
1.05, 2.39
0.029
0.10
Fair
2.02
0.86, 4.76
0.11
0.3
Poor
1.55
0.41, 5.79
0.5
0.8
Caregiver Composition
Dual-caregiver household
—
—
Single-caregiver household
1.00
0.70, 1.43
>0.9
>0.9
Neighborhood Physical Resources
0.81
0.69, 0.96
0.012
0.049
Neighborhood Environmental Threats
1.01
0.85, 1.20
>0.9
>0.9
Neighborhood Collective Efficacy
0.86
0.73, 1.02
0.086
0.2
Special Healthcare Needs
Non-SHCN
—
—
SHCN
1.78
1.25, 2.53
0.001
0.012
NDD Severity (Total)
1.15
1.08, 1.21
<0.001
<0.001
Healthcare Access Index
High Access
—
—
Medium Access
0.36
0.25, 0.51
<0.001
<0.001
Low Access
0.92
0.53, 1.59
0.8
0.9
Access Barriers Index
No Barriers
—
—
Few Barriers
0.87
0.56, 1.37
0.6
0.8
Moderate Barriers
0.89
0.38, 2.08
0.8
0.9
Many Barriers
0.75
0.19, 2.92
0.7
0.9
Insurance Type
Private Only
—
—
Private + Public
1.06
0.61, 1.84
0.8
0.9
Public Only
0.70
0.47, 1.05
0.081
0.2
Uninsured
0.09
0.01, 0.94
0.044
0.14
Insurance Adequacy Index
High Adequacy
—
—
Moderate Adequacy
1.22
0.85, 1.76
0.3
0.5
Low Adequacy
0.91
0.57, 1.46
0.7
0.9
ACEs * Child Flourishing Index
1 ACE * Flourishing
1.93
0.52, 7.23
0.3
0.6
2-3 ACEs * Flourishing
0.80
0.16, 3.92
0.8
0.9
4+ ACEs * Flourishing
1.74
0.26, 11.6
0.6
0.8
Note. This logistic regression model examines the moderating effect of the Child Flourishing Index (CFI) on the association between Adverse Childhood Experiences (ACEs) and anxiety outcomes in neurodivergent children. Survey weights from the National Survey of Children’s Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions and survey year, with FDR-adjusted p-values to account for multiple comparisons.
1OR = Odds Ratio, CI = Confidence Interval
2False discovery rate correction for multiple testing
Expand Code
# Logistic regression model with CFI_dich as a moderator for Depression, using robust standard errorsmodel_dep_mod <-svyglm(Depression ~ ace_counts * CFI_dich + sc_race_cat + sc_sex_cat + fpl_category + caregiver_education + sc_age_years + survey_year + caregiver_mental_health + child_gen_health + caregiver_composition + physical_resources + environmental_threats + collective_efficacy + special_healthcare_needs + total_severity_score + HAI_category + barriers_cat + insurance_type + IAI_category, design = design, family =quasibinomial())# Extract residuals from the logistic regression modelresiduals_dep <-residuals(model_dep_mod, type ="deviance")# Plot residuals to visually inspect themplot(residuals_dep, main ="Residuals Plot for Depression Model", ylab ="Residuals", xlab ="Index", pch =20, col ="blue")abline(h =0, col ="red", lty =2)
Expand Code
# Perform the Breusch-Pagan test for heteroscedasticitybptest(model_dep_mod)# Calculate Cook's distance for influential observations (leverage) in each modelcooks_dep <-cooks.distance(model_dep_mod)# Plot for Depression modelplot(cooks_dep, type ="p", main ="Cook's Distance for Depression Model",ylab ="Cook's Distance", xlab ="Observation Index")abline(h =4/length(cooks_dep), col ="red", lty =2) # Add a reference line for influential points
Expand Code
# Obtain robust standard errors for Depression modelrobust_dep <-coeftest(model_dep_mod, vcov =vcovHC(model_dep_mod, type ="HC0"))# Create a summary table for the Depression model with robust standard errorstbl_dep_mod <-tbl_regression(model_dep_mod, exponentiate =TRUE, label =list( ace_counts ~"ACEs", CFI_dich ~"Child Flourishing Index", fpl_category ~"Federal Poverty Level", caregiver_education ~"Caregiver Education", sc_age_years ~"Child Age", sc_race_cat ~"Child Race", sc_sex_cat ~"Child Sex", survey_year ~"Survey Year", caregiver_mental_health ~"Caregiver Mental Health", child_gen_health ~"Child General Health", caregiver_composition ~"Caregiver Composition", physical_resources ~"Neighborhood Physical Resources", environmental_threats ~"Neighborhood Environmental Threats", collective_efficacy ~"Neighborhood Collective Efficacy", special_healthcare_needs ~"Special Healthcare Needs", total_severity_score ~"NDD Severity (Total)", HAI_category ~"Healthcare Access Index", barriers_cat ~"Access Barriers Index", insurance_type ~"Insurance Type", IAI_category ~"Insurance Adequacy Index" )) %>%add_q() %>%# adjusts p-values for multiple testing (FDR)bold_p(t =0.05) %>%bold_p(t =0.05, q =TRUE) %>%# now bold q-values under the threshold of 0.05italicize_levels() %>% gtsummary::as_gt() %>% gt::tab_source_note(gt::md("*Note*. This logistic regression model examines the moderating effect of the Child Flourishing Index (CFI) on the association between Adverse Childhood Experiences (ACEs) and depression outcomes in neurodivergent children. Survey weights from the National Survey of Children's Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions and survey year, with FDR-adjusted p-values to account for multiple comparisons."))# Print table to review resultstbl_dep_mod
studentized Breusch-Pagan test
data: model_dep_mod
BP = 468.18, df = 52, p-value < 2.2e-16
Characteristic
OR1
95% CI1
p-value
q-value2
ACEs
0 ACEs
—
—
1 ACE
1.51
1.01, 2.25
0.044
0.14
2-3 ACEs
1.81
1.19, 2.74
0.005
0.039
4+ ACEs
2.15
1.39, 3.30
<0.001
0.009
Child Flourishing Index
Not Flourishing
—
—
Flourishing
1.24
0.37, 4.18
0.7
0.8
Child Race
Non-Hispanic White
—
—
Non-Hispanic Black or African American
0.76
0.46, 1.25
0.3
0.6
Non-Hispanic American Indian or Alaska Native
0.77
0.26, 2.28
0.6
0.8
Non-Hispanic Asian, Native Hawaiian, or Pacific Islander
1.29
0.63, 2.65
0.5
0.8
Non-Hispanic Other or Mixed Race
1.46
0.95, 2.24
0.081
0.2
Hispanic/Latino of Any Race
0.58
0.38, 0.88
0.010
0.058
Child Sex
Male
—
—
Female
1.43
1.07, 1.93
0.017
0.080
Federal Poverty Level
<100%
—
—
100%-199%
1.18
0.72, 1.93
0.5
0.8
200%-299%
0.80
0.44, 1.47
0.5
0.8
300%-399%
0.86
0.42, 1.75
0.7
0.8
≥400%
0.99
0.53, 1.86
>0.9
>0.9
Caregiver Education
Less Than High School
—
—
High School/Vocational
1.14
0.53, 2.41
0.7
0.8
Some College/Associate Degree
0.87
0.43, 1.79
0.7
0.8
College Degree/Higher
0.91
0.41, 2.04
0.8
0.9
Child Age
1.00
0.95, 1.05
>0.9
>0.9
Survey Year
2016
—
—
2017
1.44
0.70, 2.97
0.3
0.6
2018
1.90
1.03, 3.51
0.039
0.14
2019
1.13
0.62, 2.05
0.7
0.8
2020
1.93
1.02, 3.66
0.044
0.14
2021
1.48
0.79, 2.76
0.2
0.5
2022
1.43
0.81, 2.53
0.2
0.5
2023
1.34
0.76, 2.36
0.3
0.6
Caregiver Mental Health
Excellent
—
—
Very Good
1.05
0.74, 1.48
0.8
0.8
Good
1.94
1.30, 2.90
0.001
0.012
Fair
2.08
1.08, 3.99
0.028
0.11
Poor
3.86
1.22, 12.2
0.021
0.093
Child General Health
Excellent
—
—
Very Good
0.92
0.68, 1.25
0.6
0.8
Good
1.71
1.14, 2.57
0.009
0.058
Fair
1.57
0.78, 3.15
0.2
0.5
Poor
8.45
1.98, 36.0
0.004
0.034
Caregiver Composition
Dual-caregiver household
—
—
Single-caregiver household
0.59
0.43, 0.81
0.001
0.012
Neighborhood Physical Resources
0.93
0.78, 1.12
0.5
0.8
Neighborhood Environmental Threats
0.94
0.76, 1.15
0.5
0.8
Neighborhood Collective Efficacy
0.97
0.82, 1.15
0.7
0.8
Special Healthcare Needs
Non-SHCN
—
—
SHCN
2.16
1.48, 3.14
<0.001
0.002
NDD Severity (Total)
1.08
1.02, 1.16
0.016
0.080
Healthcare Access Index
High Access
—
—
Medium Access
0.32
0.22, 0.47
<0.001
<0.001
Low Access
0.73
0.47, 1.12
0.2
0.4
Access Barriers Index
No Barriers
—
—
Few Barriers
0.84
0.60, 1.18
0.3
0.6
Moderate Barriers
0.78
0.34, 1.76
0.5
0.8
Many Barriers
0.81
0.21, 3.07
0.8
0.8
Insurance Type
Private Only
—
—
Private + Public
1.25
0.71, 2.20
0.4
0.8
Public Only
1.07
0.71, 1.61
0.7
0.8
Uninsured
0.22
0.03, 1.57
0.13
0.4
Insurance Adequacy Index
High Adequacy
—
—
Moderate Adequacy
1.17
0.87, 1.58
0.3
0.6
Low Adequacy
1.25
0.76, 2.05
0.4
0.7
ACEs * Child Flourishing Index
1 ACE * Flourishing
0.69
0.11, 4.38
0.7
0.8
2-3 ACEs * Flourishing
4.35
0.88, 21.4
0.070
0.2
4+ ACEs * Flourishing
1.50
0.30, 7.44
0.6
0.8
Note. This logistic regression model examines the moderating effect of the Child Flourishing Index (CFI) on the association between Adverse Childhood Experiences (ACEs) and depression outcomes in neurodivergent children. Survey weights from the National Survey of Children’s Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions and survey year, with FDR-adjusted p-values to account for multiple comparisons.
1OR = Odds Ratio, CI = Confidence Interval
2False discovery rate correction for multiple testing
Expand Code
# Logistic regression model with CFI_dich as a moderator for Behavioral Problems, using robust standard errorsmodel_beh_mod <-svyglm(Behavioral_Problems ~ ace_counts * CFI_dich + sc_race_cat + sc_sex_cat + fpl_category + caregiver_education + sc_age_years + survey_year + caregiver_mental_health + child_gen_health + caregiver_composition + physical_resources + environmental_threats + collective_efficacy + special_healthcare_needs + total_severity_score + HAI_category + barriers_cat + insurance_type + IAI_category,design = design, family =quasibinomial())# Extract residuals from the logistic regression modelresiduals_beh <-residuals(model_beh_mod, type ="deviance")# Plot residuals to visually inspect themplot(residuals_beh, main ="Residuals Plot for Behavioral Problems Model", ylab ="Residuals", xlab ="Index", pch =20, col ="blue")abline(h =0, col ="red", lty =2)
Expand Code
# Perform the Breusch-Pagan test for heteroscedasticitybptest(model_beh_mod)# Calculate Cook's distance for influential observations (leverage) in each modelcooks_beh <-cooks.distance(model_beh_mod)# Plot for Behavioral Problems modelplot(cooks_beh, type ="p", main ="Cook's Distance for Behavioral Problems Model",ylab ="Cook's Distance", xlab ="Observation Index")abline(h =4/length(cooks_beh), col ="red", lty =2) # Add a reference line for influential points
Expand Code
# Obtain robust standard errors for Behavioral Problems modelrobust_beh <-coeftest(model_beh_mod, vcov =vcovHC(model_beh_mod, type ="HC0"))# Create a summary table for the Behavioral Problems model with robust standard errorstbl_beh_mod <-tbl_regression(model_beh_mod, exponentiate =TRUE, label =list( ace_counts ~"ACEs", CFI_dich ~"Child Flourishing Index", fpl_category ~"Federal Poverty Level", caregiver_education ~"Caregiver Education", sc_age_years ~"Child Age", sc_race_cat ~"Child Race", sc_sex_cat ~"Child Sex", survey_year ~"Survey Year", caregiver_mental_health ~"Caregiver Mental Health", child_gen_health ~"Child General Health", caregiver_composition ~"Caregiver Composition", physical_resources ~"Neighborhood Physical Resources", environmental_threats ~"Neighborhood Environmental Threats", collective_efficacy ~"Neighborhood Collective Efficacy", special_healthcare_needs ~"Special Healthcare Needs", total_severity_score ~"NDD Severity (Total)", HAI_category ~"Healthcare Access Index", barriers_cat ~"Access Barriers Index", insurance_type ~"Insurance Type", IAI_category ~"Insurance Adequacy Index" )) %>%add_q() %>%# adjusts p-values for multiple testing (FDR)bold_p(t =0.05) %>%bold_p(t =0.05, q =TRUE) %>%# now bold q-values under the threshold of 0.05italicize_levels() %>% gtsummary::as_gt() %>% gt::tab_source_note(gt::md("*Note*. This logistic regression model examines the moderating effect of the Child Flourishing Index (CFI) on the association between Adverse Childhood Experiences (ACEs) and behavioral problem outcomes in neurodivergent children. Survey weights from the National Survey of Children's Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions and survey year, with FDR-adjusted p-values to account for multiple comparisons."))# Print table to review resultstbl_beh_mod
studentized Breusch-Pagan test
data: model_beh_mod
BP = 1219.8, df = 52, p-value < 2.2e-16
Characteristic
OR1
95% CI1
p-value
q-value2
ACEs
0 ACEs
—
—
1 ACE
1.04
0.77, 1.41
0.8
>0.9
2-3 ACEs
0.84
0.61, 1.18
0.3
0.8
4+ ACEs
1.19
0.82, 1.73
0.4
0.8
Child Flourishing Index
Not Flourishing
—
—
Flourishing
0.79
0.36, 1.72
0.5
0.9
Child Race
Non-Hispanic White
—
—
Non-Hispanic Black or African American
1.06
0.75, 1.48
0.8
0.9
Non-Hispanic American Indian or Alaska Native
0.80
0.33, 1.92
0.6
0.9
Non-Hispanic Asian, Native Hawaiian, or Pacific Islander
0.61
0.31, 1.21
0.2
0.5
Non-Hispanic Other or Mixed Race
0.93
0.67, 1.29
0.7
0.9
Hispanic/Latino of Any Race
0.90
0.67, 1.22
0.5
0.9
Child Sex
Male
—
—
Female
1.04
0.83, 1.30
0.8
0.9
Federal Poverty Level
<100%
—
—
100%-199%
0.92
0.63, 1.32
0.6
0.9
200%-299%
1.12
0.74, 1.69
0.6
0.9
300%-399%
0.81
0.52, 1.27
0.4
0.8
≥400%
0.90
0.58, 1.39
0.6
0.9
Caregiver Education
Less Than High School
—
—
High School/Vocational
0.83
0.42, 1.63
0.6
0.9
Some College/Associate Degree
0.77
0.39, 1.50
0.4
0.8
College Degree/Higher
0.67
0.34, 1.32
0.2
0.7
Child Age
0.83
0.80, 0.86
<0.001
<0.001
Survey Year
2016
—
—
2017
0.80
0.53, 1.22
0.3
0.8
2018
0.78
0.54, 1.13
0.2
0.6
2019
0.92
0.61, 1.41
0.7
0.9
2020
1.40
0.91, 2.14
0.12
0.5
2021
1.02
0.70, 1.48
>0.9
>0.9
2022
1.07
0.75, 1.52
0.7
0.9
2023
0.74
0.52, 1.06
0.11
0.5
Caregiver Mental Health
Excellent
—
—
Very Good
1.08
0.84, 1.39
0.5
0.9
Good
1.00
0.74, 1.35
>0.9
>0.9
Fair
1.12
0.69, 1.80
0.6
0.9
Poor
1.15
0.34, 3.92
0.8
>0.9
Child General Health
Excellent
—
—
Very Good
1.33
1.05, 1.67
0.018
0.12
Good
1.29
0.96, 1.74
0.10
0.5
Fair
2.49
1.40, 4.40
0.002
0.015
Poor
26.7
5.54, 129
<0.001
<0.001
Caregiver Composition
Dual-caregiver household
—
—
Single-caregiver household
0.90
0.70, 1.14
0.4
0.8
Neighborhood Physical Resources
1.00
0.90, 1.12
>0.9
>0.9
Neighborhood Environmental Threats
0.96
0.85, 1.09
0.5
0.9
Neighborhood Collective Efficacy
1.11
0.99, 1.25
0.084
0.4
Special Healthcare Needs
Non-SHCN
—
—
SHCN
3.23
2.46, 4.23
<0.001
<0.001
NDD Severity (Total)
1.17
1.12, 1.22
<0.001
<0.001
Healthcare Access Index
High Access
—
—
Medium Access
0.67
0.53, 0.84
<0.001
0.007
Low Access
1.03
0.72, 1.48
0.9
>0.9
Access Barriers Index
No Barriers
—
—
Few Barriers
0.94
0.66, 1.33
0.7
0.9
Moderate Barriers
0.74
0.37, 1.51
0.4
0.8
Many Barriers
1.07
0.29, 3.90
>0.9
>0.9
Insurance Type
Private Only
—
—
Private + Public
1.38
0.89, 2.15
0.2
0.5
Public Only
1.07
0.81, 1.41
0.6
0.9
Uninsured
9.45
0.94, 95.0
0.057
0.3
Insurance Adequacy Index
High Adequacy
—
—
Moderate Adequacy
0.70
0.55, 0.90
0.005
0.037
Low Adequacy
0.86
0.62, 1.18
0.4
0.8
ACEs * Child Flourishing Index
1 ACE * Flourishing
0.50
0.11, 2.20
0.4
0.8
2-3 ACEs * Flourishing
2.37
0.61, 9.22
0.2
0.7
4+ ACEs * Flourishing
0.51
0.13, 1.99
0.3
0.8
Note. This logistic regression model examines the moderating effect of the Child Flourishing Index (CFI) on the association between Adverse Childhood Experiences (ACEs) and behavioral problem outcomes in neurodivergent children. Survey weights from the National Survey of Children’s Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions and survey year, with FDR-adjusted p-values to account for multiple comparisons.
1OR = Odds Ratio, CI = Confidence Interval
2False discovery rate correction for multiple testing
The anxiety model did not yield significant interactions, so we do not proceed with stratified post-hoc analyses. Instead, we will run main effect models without interaction terms.
The depression model did not yield significant interactions, so we do not proceed with stratified post-hoc analyses. Instead, we will run main effect models without interaction terms.
The behavioral problems model didnot yield significant interactions, so we do not proceed stratified with post-hoc analyses. Instead, we will run main effect models without interaction terms.
Main Effect Models
Fit three logistic regression models with survey weights to examine the main effects of ACE counts and the Child Flourishing Index (CFI_dich) on psychopathology outcomes (anxiety, depression, behavioral problems). Extract residuals, check heteroscedasticity using Breusch-Pagan tests, calculate Cook’s distance for influential observations, and obtain robust standard errors. Create summary tables with robust standard errors and False Discovery Rate (FDR) adjusted p-values for each model.
# Logistic regression model for main effects using robust standard errorsmodel_anx <-svyglm(Anxiety ~ ace_counts + CFI_dich + sc_race_cat + sc_sex_cat + fpl_category + caregiver_education + sc_age_years + survey_year + caregiver_mental_health + child_gen_health + caregiver_composition + physical_resources + environmental_threats + collective_efficacy + special_healthcare_needs + total_severity_score + HAI_category + barriers_cat + insurance_type + IAI_category, design = design, family =quasibinomial())# Extract residuals from the logistic regression modelresiduals_anx <-residuals(model_anx, type ="deviance")# Plot residuals to visually inspect themplot(residuals_anx, main ="Residuals Plot for Anxiety Model", ylab ="Residuals", xlab ="Index", pch =20, col ="blue")abline(h =0, col ="red", lty =2)
Expand Code
# Perform the Breusch-Pagan test for heteroscedasticitybptest(model_anx)
studentized Breusch-Pagan test
data: model_anx
BP = 849.67, df = 49, p-value < 2.2e-16
Expand Code
# Calculate Cook's distance for influential observations (leverage) in each modelcooks_anx <-cooks.distance(model_anx)# Plot for Anxiety modelplot(cooks_anx, type ="p", main ="Cook's Distance for Anxiety Model",ylab ="Cook's Distance", xlab ="Observation Index")abline(h =4/length(cooks_anx), col ="red", lty =2) # Add a reference line for influential points
Expand Code
# Obtain robust standard errors for Anxiety modelrobust_anx <-coeftest(model_anx, vcov =vcovHC(model_anx, type ="HC0"))# Create a summary table for the Anxiety model with robust standard errorstbl_anx <-tbl_regression(model_anx, exponentiate =TRUE, label =list( ace_counts ~"ACEs", CFI_dich ~"Child Flourishing Index", fpl_category ~"Federal Poverty Level", caregiver_education ~"Caregiver Education", sc_age_years ~"Child Age", sc_race_cat ~"Child Race", sc_sex_cat ~"Child Sex", survey_year ~"Survey Year", caregiver_mental_health ~"Caregiver Mental Health", child_gen_health ~"Child General Health", caregiver_composition ~"Caregiver Composition", physical_resources ~"Neighborhood Physical Resources", environmental_threats ~"Neighborhood Environmental Threats", collective_efficacy ~"Neighborhood Collective Efficacy", special_healthcare_needs ~"Special Healthcare Needs", total_severity_score ~"NDD Severity (Total)", HAI_category ~"Healthcare Access Index", barriers_cat ~"Access Barriers Index", insurance_type ~"Insurance Type", IAI_category ~"Insurance Adequacy Index" )) %>%add_q() %>%# adjusts p-values for multiple testing (FDR)bold_p(t =0.05) %>%bold_p(t =0.05, q =TRUE) %>%# now bold q-values under the threshold of 0.05italicize_levels() %>% gtsummary::as_gt() %>% gt::tab_source_note(gt::md("*Note*. This logistic regression model examines the association between Adverse Childhood Experiences (ACEs) and anxiety outcomes in neurodivergent children, including main effects of the Child Flourishing Index (CFI) without interaction terms. Survey weights from the National Survey of Children's Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions, and survey year, with FDR-adjusted p-values to account for multiple comparisons."))# Print table to review resultstbl_anx
Characteristic
OR1
95% CI1
p-value
q-value2
ACEs
0 ACEs
—
—
1 ACE
1.15
0.77, 1.72
0.5
0.7
2-3 ACEs
1.38
0.89, 2.14
0.2
0.3
4+ ACEs
1.56
0.96, 2.54
0.074
0.2
Child Flourishing Index
Not Flourishing
—
—
Flourishing
1.27
0.62, 2.63
0.5
0.8
Child Race
Non-Hispanic White
—
—
Non-Hispanic Black or African American
0.78
0.45, 1.36
0.4
0.6
Non-Hispanic American Indian or Alaska Native
0.58
0.23, 1.44
0.2
0.4
Non-Hispanic Asian, Native Hawaiian, or Pacific Islander
0.33
0.17, 0.65
0.001
0.010
Non-Hispanic Other or Mixed Race
0.95
0.61, 1.46
0.8
0.9
Hispanic/Latino of Any Race
0.61
0.43, 0.86
0.005
0.021
Child Sex
Male
—
—
Female
0.98
0.73, 1.32
0.9
>0.9
Federal Poverty Level
<100%
—
—
100%-199%
0.59
0.35, 1.00
0.050
0.14
200%-299%
0.77
0.43, 1.37
0.4
0.6
300%-399%
0.51
0.28, 0.93
0.028
0.093
≥400%
0.75
0.42, 1.34
0.3
0.6
Caregiver Education
Less Than High School
—
—
High School/Vocational
3.09
1.41, 6.75
0.005
0.021
Some College/Associate Degree
2.97
1.43, 6.16
0.003
0.021
College Degree/Higher
4.20
1.98, 8.92
<0.001
0.003
Child Age
0.99
0.95, 1.04
0.7
0.9
Survey Year
2016
—
—
2017
1.12
0.60, 2.09
0.7
0.9
2018
0.92
0.56, 1.52
0.7
0.9
2019
0.64
0.39, 1.06
0.084
0.2
2020
1.20
0.67, 2.14
0.5
0.8
2021
1.90
1.22, 2.96
0.005
0.021
2022
1.44
0.89, 2.33
0.14
0.3
2023
1.23
0.76, 2.01
0.4
0.6
Caregiver Mental Health
Excellent
—
—
Very Good
1.24
0.90, 1.72
0.2
0.4
Good
1.76
1.21, 2.58
0.003
0.021
Fair
2.43
1.42, 4.17
0.001
0.010
Poor
5.80
1.72, 19.6
0.005
0.021
Child General Health
Excellent
—
—
Very Good
1.04
0.77, 1.39
0.8
0.9
Good
1.59
1.05, 2.40
0.027
0.093
Fair
2.02
0.86, 4.76
0.11
0.2
Poor
1.52
0.39, 5.88
0.5
0.8
Caregiver Composition
Dual-caregiver household
—
—
Single-caregiver household
1.00
0.70, 1.43
>0.9
>0.9
Neighborhood Physical Resources
0.81
0.69, 0.95
0.011
0.042
Neighborhood Environmental Threats
1.01
0.85, 1.20
>0.9
>0.9
Neighborhood Collective Efficacy
0.87
0.74, 1.02
0.092
0.2
Special Healthcare Needs
Non-SHCN
—
—
SHCN
1.78
1.26, 2.53
0.001
0.010
NDD Severity (Total)
1.14
1.08, 1.21
<0.001
<0.001
Healthcare Access Index
High Access
—
—
Medium Access
0.36
0.25, 0.51
<0.001
<0.001
Low Access
0.92
0.53, 1.60
0.8
0.9
Access Barriers Index
No Barriers
—
—
Few Barriers
0.88
0.56, 1.37
0.6
0.8
Moderate Barriers
0.90
0.38, 2.08
0.8
0.9
Many Barriers
0.76
0.20, 2.97
0.7
0.9
Insurance Type
Private Only
—
—
Private + Public
1.06
0.61, 1.84
0.8
0.9
Public Only
0.70
0.46, 1.04
0.078
0.2
Uninsured
0.09
0.01, 0.93
0.044
0.13
Insurance Adequacy Index
High Adequacy
—
—
Moderate Adequacy
1.22
0.85, 1.75
0.3
0.5
Low Adequacy
0.91
0.57, 1.46
0.7
0.9
Note. This logistic regression model examines the association between Adverse Childhood Experiences (ACEs) and anxiety outcomes in neurodivergent children, including main effects of the Child Flourishing Index (CFI) without interaction terms. Survey weights from the National Survey of Children’s Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions, and survey year, with FDR-adjusted p-values to account for multiple comparisons.
1OR = Odds Ratio, CI = Confidence Interval
2False discovery rate correction for multiple testing
Expand Code
# Logistic regression model for main effects using robust standard errorsmodel_dep <-svyglm(Depression ~ ace_counts + CFI_dich + sc_race_cat + sc_sex_cat + fpl_category + caregiver_education + sc_age_years + survey_year + caregiver_mental_health + child_gen_health + caregiver_composition + physical_resources + environmental_threats + collective_efficacy + special_healthcare_needs + total_severity_score + HAI_category + barriers_cat + insurance_type + IAI_category, design = design, family =quasibinomial())# Extract residuals from the logistic regression modelresiduals_dep <-residuals(model_dep, type ="deviance")# Plot residuals to visually inspect themplot(residuals_dep, main ="Residuals Plot for Depression Model", ylab ="Residuals", xlab ="Index", pch =20, col ="blue")abline(h =0, col ="red", lty =2)
Expand Code
# Perform the Breusch-Pagan test for heteroscedasticitybptest(model_dep)
studentized Breusch-Pagan test
data: model_dep
BP = 467.39, df = 49, p-value < 2.2e-16
Expand Code
# Calculate Cook's distance for influential observations (leverage) in each modelcooks_dep <-cooks.distance(model_dep)# Plot for Depression modelplot(cooks_dep, type ="p", main ="Cook's Distance for Depression Model",ylab ="Cook's Distance", xlab ="Observation Index")abline(h =4/length(cooks_dep), col ="red", lty =2) # Add a reference line for influential points
Expand Code
# Obtain robust standard errors for Depression modelrobust_dep <-coeftest(model_dep, vcov =vcovHC(model_dep, type ="HC0"))# Create a summary table for the Depression model with robust standard errorstbl_dep <-tbl_regression(model_dep, exponentiate =TRUE, label =list( ace_counts ~"ACEs", CFI_dich ~"Child Flourishing Index", fpl_category ~"Federal Poverty Level", caregiver_education ~"Caregiver Education", sc_age_years ~"Child Age", sc_race_cat ~"Child Race", sc_sex_cat ~"Child Sex", survey_year ~"Survey Year", caregiver_mental_health ~"Caregiver Mental Health", child_gen_health ~"Child General Health", caregiver_composition ~"Caregiver Composition", physical_resources ~"Neighborhood Physical Resources", environmental_threats ~"Neighborhood Environmental Threats", collective_efficacy ~"Neighborhood Collective Efficacy", special_healthcare_needs ~"Special Healthcare Needs", total_severity_score ~"NDD Severity (Total)", HAI_category ~"Healthcare Access Index", barriers_cat ~"Access Barriers Index", insurance_type ~"Insurance Type", IAI_category ~"Insurance Adequacy Index" )) %>%add_q() %>%# adjusts p-values for multiple testing (FDR)bold_p(t =0.05) %>%bold_p(t =0.05, q =TRUE) %>%# now bold q-values under the threshold of 0.05italicize_levels() %>% gtsummary::as_gt() %>% gt::tab_source_note(gt::md("*Note*. This logistic regression model examines the association between Adverse Childhood Experiences (ACEs) and depression outcomes in neurodivergent children, including main effects of the Child Flourishing Index (CFI) without interaction terms. Survey weights from the National Survey of Children's Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions, and survey year, with FDR-adjusted p-values to account for multiple comparisons."))# Print table to review resultstbl_dep
Characteristic
OR1
95% CI1
p-value
q-value2
ACEs
0 ACEs
—
—
1 ACE
1.50
1.02, 2.23
0.042
0.13
2-3 ACEs
1.85
1.23, 2.79
0.003
0.027
4+ ACEs
2.17
1.41, 3.33
<0.001
0.006
Child Flourishing Index
Not Flourishing
—
—
Flourishing
1.84
0.96, 3.52
0.065
0.2
Child Race
Non-Hispanic White
—
—
Non-Hispanic Black or African American
0.76
0.46, 1.24
0.3
0.5
Non-Hispanic American Indian or Alaska Native
0.77
0.26, 2.29
0.6
0.8
Non-Hispanic Asian, Native Hawaiian, or Pacific Islander
1.29
0.63, 2.65
0.5
0.7
Non-Hispanic Other or Mixed Race
1.45
0.95, 2.22
0.085
0.2
Hispanic/Latino of Any Race
0.58
0.39, 0.88
0.010
0.056
Child Sex
Male
—
—
Female
1.42
1.06, 1.91
0.019
0.083
Federal Poverty Level
<100%
—
—
100%-199%
1.17
0.71, 1.91
0.5
0.7
200%-299%
0.79
0.43, 1.46
0.5
0.7
300%-399%
0.85
0.42, 1.73
0.7
0.8
≥400%
0.98
0.52, 1.84
>0.9
>0.9
Caregiver Education
Less Than High School
—
—
High School/Vocational
1.15
0.55, 2.44
0.7
0.8
Some College/Associate Degree
0.89
0.44, 1.81
0.7
0.8
College Degree/Higher
0.93
0.42, 2.07
0.9
0.9
Child Age
1.00
0.95, 1.05
>0.9
>0.9
Survey Year
2016
—
—
2017
1.46
0.71, 3.00
0.3
0.5
2018
1.91
1.04, 3.51
0.037
0.13
2019
1.14
0.63, 2.07
0.7
0.8
2020
1.94
1.02, 3.67
0.042
0.13
2021
1.49
0.80, 2.78
0.2
0.5
2022
1.45
0.82, 2.56
0.2
0.5
2023
1.35
0.77, 2.38
0.3
0.5
Caregiver Mental Health
Excellent
—
—
Very Good
1.05
0.74, 1.48
0.8
0.8
Good
1.95
1.31, 2.90
0.001
0.010
Fair
2.08
1.09, 3.98
0.027
0.10
Poor
3.89
1.23, 12.3
0.021
0.085
Child General Health
Excellent
—
—
Very Good
0.92
0.68, 1.25
0.6
0.8
Good
1.71
1.14, 2.56
0.010
0.056
Fair
1.57
0.78, 3.16
0.2
0.5
Poor
8.48
1.99, 36.2
0.004
0.027
Caregiver Composition
Dual-caregiver household
—
—
Single-caregiver household
0.59
0.43, 0.81
<0.001
0.010
Neighborhood Physical Resources
0.93
0.78, 1.12
0.5
0.7
Neighborhood Environmental Threats
0.94
0.76, 1.16
0.6
0.7
Neighborhood Collective Efficacy
0.97
0.82, 1.15
0.7
0.8
Special Healthcare Needs
Non-SHCN
—
—
SHCN
2.14
1.47, 3.12
<0.001
0.002
NDD Severity (Total)
1.08
1.02, 1.16
0.016
0.078
Healthcare Access Index
High Access
—
—
Medium Access
0.32
0.22, 0.46
<0.001
<0.001
Low Access
0.73
0.47, 1.12
0.2
0.4
Access Barriers Index
No Barriers
—
—
Few Barriers
0.84
0.60, 1.18
0.3
0.5
Moderate Barriers
0.78
0.34, 1.76
0.5
0.7
Many Barriers
0.78
0.21, 2.93
0.7
0.8
Insurance Type
Private Only
—
—
Private + Public
1.23
0.70, 2.17
0.5
0.7
Public Only
1.07
0.71, 1.61
0.7
0.8
Uninsured
0.23
0.03, 1.58
0.13
0.3
Insurance Adequacy Index
High Adequacy
—
—
Moderate Adequacy
1.17
0.87, 1.58
0.3
0.5
Low Adequacy
1.25
0.76, 2.04
0.4
0.6
Note. This logistic regression model examines the association between Adverse Childhood Experiences (ACEs) and depression outcomes in neurodivergent children, including main effects of the Child Flourishing Index (CFI) without interaction terms. Survey weights from the National Survey of Children’s Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions, and survey year, with FDR-adjusted p-values to account for multiple comparisons.
1OR = Odds Ratio, CI = Confidence Interval
2False discovery rate correction for multiple testing
Expand Code
# Logistic regression model for main effects using robust standard errorsmodel_beh <-svyglm(Behavioral_Problems ~ ace_counts + CFI_dich + sc_race_cat + sc_sex_cat + fpl_category + caregiver_education + sc_age_years + survey_year + caregiver_mental_health + child_gen_health + caregiver_composition + physical_resources + environmental_threats + collective_efficacy + special_healthcare_needs + total_severity_score + HAI_category + barriers_cat + insurance_type + IAI_category,design = design, family =quasibinomial())# Extract residuals from the logistic regression modelresiduals_beh <-residuals(model_beh, type ="deviance")# Plot residuals to visually inspect themplot(residuals_beh, main ="Residuals Plot for Behavioral Problems Model", ylab ="Residuals", xlab ="Index", pch =20, col ="blue")abline(h =0, col ="red", lty =2)
Expand Code
# Perform the Breusch-Pagan test for heteroscedasticitybptest(model_beh)
studentized Breusch-Pagan test
data: model_beh
BP = 1228.6, df = 49, p-value < 2.2e-16
Expand Code
# Calculate Cook's distance for influential observations (leverage) in each modelcooks_beh <-cooks.distance(model_beh)# Plot for Behavioral Problems modelplot(cooks_beh, type ="p", main ="Cook's Distance for Behavioral Problems Model",ylab ="Cook's Distance", xlab ="Observation Index")abline(h =4/length(cooks_beh), col ="red", lty =2) # Add a reference line for influential points
Expand Code
# Obtain robust standard errors for Behavioral Problems modelrobust_beh <-coeftest(model_beh, vcov =vcovHC(model_beh, type ="HC0"))# Create a summary table for the Behavioral Problems model with robust standard errorstbl_beh <-tbl_regression(model_beh, exponentiate =TRUE, label =list( ace_counts ~"ACEs", CFI_dich ~"Child Flourishing Index", fpl_category ~"Federal Poverty Level", caregiver_education ~"Caregiver Education", sc_age_years ~"Child Age", sc_race_cat ~"Child Race", sc_sex_cat ~"Child Sex", survey_year ~"Survey Year", caregiver_mental_health ~"Caregiver Mental Health", child_gen_health ~"Child General Health", caregiver_composition ~"Caregiver Composition", physical_resources ~"Neighborhood Physical Resources", environmental_threats ~"Neighborhood Environmental Threats", collective_efficacy ~"Neighborhood Collective Efficacy", special_healthcare_needs ~"Special Healthcare Needs", total_severity_score ~"NDD Severity (Total)", HAI_category ~"Healthcare Access Index", barriers_cat ~"Access Barriers Index", insurance_type ~"Insurance Type", IAI_category ~"Insurance Adequacy Index" )) %>%add_q() %>%# adjusts p-values for multiple testing (FDR)bold_p(t =0.05) %>%bold_p(t =0.05, q =TRUE) %>%# now bold q-values under the threshold of 0.05italicize_levels() %>% gtsummary::as_gt() %>% gt::tab_source_note(gt::md("*Note*. This logistic regression model examines the association between Adverse Childhood Experiences (ACEs) and behavioral problem outcomes in neurodivergent children, including main effects of the Child Flourishing Index (CFI) without interaction terms. Survey weights from the National Survey of Children's Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions, and survey year, with FDR-adjusted p-values to account for multiple comparisons."))# Print table to review resultstbl_beh
Characteristic
OR1
95% CI1
p-value
q-value2
ACEs
0 ACEs
—
—
1 ACE
1.00
0.74, 1.37
>0.9
>0.9
2-3 ACEs
0.86
0.62, 1.19
0.4
0.8
4+ ACEs
1.18
0.81, 1.71
0.4
0.8
Child Flourishing Index
Not Flourishing
—
—
Flourishing
0.66
0.27, 1.59
0.4
0.8
Child Race
Non-Hispanic White
—
—
Non-Hispanic Black or African American
1.06
0.76, 1.49
0.7
0.9
Non-Hispanic American Indian or Alaska Native
0.81
0.33, 1.95
0.6
0.9
Non-Hispanic Asian, Native Hawaiian, or Pacific Islander
0.61
0.31, 1.20
0.2
0.5
Non-Hispanic Other or Mixed Race
0.95
0.68, 1.33
0.8
0.9
Hispanic/Latino of Any Race
0.90
0.66, 1.22
0.5
0.9
Child Sex
Male
—
—
Female
1.04
0.83, 1.31
0.7
0.9
Federal Poverty Level
<100%
—
—
100%-199%
0.91
0.63, 1.31
0.6
0.9
200%-299%
1.11
0.74, 1.68
0.6
0.9
300%-399%
0.80
0.51, 1.25
0.3
0.8
≥400%
0.88
0.57, 1.36
0.6
0.9
Caregiver Education
Less Than High School
—
—
High School/Vocational
0.83
0.43, 1.63
0.6
0.9
Some College/Associate Degree
0.77
0.40, 1.51
0.5
0.9
College Degree/Higher
0.68
0.34, 1.35
0.3
0.8
Child Age
0.83
0.80, 0.86
<0.001
<0.001
Survey Year
2016
—
—
2017
0.81
0.53, 1.24
0.3
0.8
2018
0.79
0.54, 1.14
0.2
0.6
2019
0.92
0.60, 1.40
0.7
0.9
2020
1.38
0.91, 2.11
0.13
0.5
2021
1.01
0.70, 1.47
>0.9
>0.9
2022
1.06
0.74, 1.50
0.8
0.9
2023
0.74
0.52, 1.06
0.10
0.5
Caregiver Mental Health
Excellent
—
—
Very Good
1.08
0.83, 1.40
0.6
0.9
Good
0.99
0.73, 1.34
>0.9
>0.9
Fair
1.11
0.69, 1.79
0.7
0.9
Poor
1.14
0.33, 3.88
0.8
>0.9
Child General Health
Excellent
—
—
Very Good
1.33
1.05, 1.69
0.016
0.10
Good
1.30
0.96, 1.75
0.090
0.4
Fair
2.51
1.41, 4.44
0.002
0.013
Poor
26.7
5.54, 128
<0.001
<0.001
Caregiver Composition
Dual-caregiver household
—
—
Single-caregiver household
0.90
0.70, 1.15
0.4
0.8
Neighborhood Physical Resources
1.00
0.90, 1.11
>0.9
>0.9
Neighborhood Environmental Threats
0.96
0.85, 1.09
0.5
0.9
Neighborhood Collective Efficacy
1.10
0.98, 1.24
0.12
0.5
Special Healthcare Needs
Non-SHCN
—
—
SHCN
3.20
2.44, 4.20
<0.001
<0.001
NDD Severity (Total)
1.17
1.12, 1.22
<0.001
<0.001
Healthcare Access Index
High Access
—
—
Medium Access
0.68
0.54, 0.85
<0.001
0.009
Low Access
1.04
0.73, 1.49
0.8
>0.9
Access Barriers Index
No Barriers
—
—
Few Barriers
0.94
0.66, 1.33
0.7
0.9
Moderate Barriers
0.74
0.37, 1.50
0.4
0.8
Many Barriers
1.06
0.29, 3.88
>0.9
>0.9
Insurance Type
Private Only
—
—
Private + Public
1.38
0.89, 2.16
0.2
0.5
Public Only
1.06
0.81, 1.40
0.7
0.9
Uninsured
9.49
0.95, 95.2
0.056
0.3
Insurance Adequacy Index
High Adequacy
—
—
Moderate Adequacy
0.70
0.54, 0.89
0.005
0.032
Low Adequacy
0.85
0.62, 1.17
0.3
0.8
Note. This logistic regression model examines the association between Adverse Childhood Experiences (ACEs) and behavioral problem outcomes in neurodivergent children, including main effects of the Child Flourishing Index (CFI) without interaction terms. Survey weights from the National Survey of Children’s Health (NSCH) were applied to ensure the representativeness of the target population. Robust standard errors (HC0) were used to account for potential heteroscedasticity not fully addressed by the survey weights. HC0 provides the original White standard errors, which are suitable for large samples and do not include small-sample bias corrections. Covariates include child demographics, family characteristics, neighborhood conditions, and survey year, with FDR-adjusted p-values to account for multiple comparisons.
1OR = Odds Ratio, CI = Confidence Interval
2False discovery rate correction for multiple testing
Session Information
To enhance reproducibility, details about the working environment used for these analyses can be found below.